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ABSTRACT 

Steady state distribution of the probabihty amphtudes and the site popula- 
tion in the one dimensional optical superlattice was found. It was shown that 
this solution of the equations which describe the dynamics of the Bose- Einstein 
condensate in superlattice is unstable at the sufficiently high density of the 
bosons. The expression for increment of the modulational instability was found 
on base of the linear stability analysis. Numerical simulation demonstrates the 
evolution of the steady state distributions of bosons into the space array of the 
solitary peaks before the chaotic regime generation. 
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1 Introduction 



Interference of the several plane waves of monochromatic radiation can form 
a diffraction pattern in which the electric field strength periodically varies in 
space. Resulting periodic system of the microscopic potentials is designated as 
optical lattice [T1I21I3]- Bose-Einstein condensates (BECs) trapped in optical 
lattice have been studied in sufficient detail [H O [5] . The dynamics of atoms in 
an optical lattice can be described by two methods. The first method is based 
on the nonlinear Schrodinger equation with a periodic potential [71 [H] • In papers 
devoted to investigation of the Bose-Einstein condensation, this equation is often 
called the Gross-Pitaevskii equation [21 UHl HI] ■ The second method answering 
the tight-binding approximation is based on the Hubbard model, where the 
operators of creation and annihilation of fermions are replaced by the operators 
satisfying the commutation relations for bosons. The resulting model is often 
called the Bose-Hubbard model [H O O [15] . 

Recently the optical lattices with two sorts of microscopic potentials in a 
unit cell have attracted attention [TBI HZl Hll HH HH] • In analogy with theory 
of solid state crystals, this periodic system of the microscopic potentials can 
be named optical superlattice. Due to the difference of the energy levels of 
microscopic potentials quantum tunneling between the nearest-neighbor sites is 
absent. Also, at low temperatures thermoactivated transport of an atom from 
one site of the micropotential to the other site is absent. However, the photo- 
induced transport of atoms along the optical superlattice is possible under the 
condition for the Raman resonance [12[ I21[ I22j . A similar process is known in 
nonlinear optics as a coherent population transfer and in solid state physics of 
low dimensional systems as coherent transfer of electrons or excitons in a system 
of coupled quantum dots. 

Under CW electromagnetic radiation the different nonlinear excitations can 
propagate in the optical superlattice. It is important to emphasize that parame- 
ters of these excitations can be controlled by additional radiation, which defines 
hopping rates between the adjacent sites. 
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Frequently the steady state wave motion takes place in various physical 
systems. Due to the interplay between dispersive and nonlincarity effects a 
week perturbation of the steady state wave may induce the exponential growth 
of the perturbation. That phenomenon is the modulational instability of a 
steady state wave. Some times modulational instability results in train of the 
soliton-like waves. But it does not always happen. 

Modulational instability in BECs in the case of the ordinary optical lattices 
was investigated in [231 IMl US . It was shown that modulational instability is 
basic mechanism by which solitons are created in BEC. The kind of solitons 
(i.e., bright or dark one) depends on sign of the scattering lengths. 

In this paper we consider the dynamics of the site populations for optical 
superlattice with two sorts of microscopic potentials in a unit cell. The generated 
Bose-Hubbard model describing this system in the tight-binding approximation 
was used in |211 122] to write the system of equations of motion for probability 
amplitudes of the site population. These equations are employed as the basis for 
present investigation. We found the stationary distribution of the probability 
amplitudes. It should be pointed out that this distribution is not homogeneous 
one, as opposed to the case of ordinary optical lattices. The main result is an 
analytic expression for the modulational instability increment, which depends on 
the site population and the wave number of harmonic weak perturbations. This 
expression was obtained on base of the linear stability analysis. The numerical 
simulation shows that instability leads to strongly non-regular pattern. 

2 Model and basic equations 

Let us consider the one-dimensional optical superlattice with two kinds of sites 
in the tight-binding approximation [121 113) . The microscopic potentials of one 
depth correspond to sites with an even number, and microscopic potentials of 
other depth correspond to sites with an odd number. Sites labeled by even 
numbers contain bosons in the ground state \ga > with the energy Sa- Sites 
labeled by odd numbers contain bosons in the ground state \gb > with the 
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energy £& > Sa- Let the temperature of the system be such that the higher 
levels of the microscopic potentials are not populated. Since the energies of 
the ground states of neighboring sites differ, the process of direct tunneling of 
an atom from one site to another one can be excluded from the consideration. 
Let us assume that biharmonic radiation (wi and L02 are the frequencies of the 
carrier waves) acts on the atoms and the condition for the Raman resonance 
(sf, — £a) ~ h{uji — UJ2) is fulfilled. In this case, after absorption of the first 
(second) photon, the atom will go from the deep (shallow) microscopic potential 
to the state of the continuous spectrum, and, after emission of the second (first) 
photon, it will return to the state having the energy Sb (ca) and thus will be 
brought into a shallow (deep) microscopic potential. Thus, although tunneling 
or thcrmoactivated transport of atoms along the optical superlattice is absent, 
their photo-induced transport is possible. 

In classical limit the system of equations describing the probability ampli- 
tudes of the populations for sites a2j — (a2j) and fc2j"+i = (^2^+1) takes the 
following form [HI [22]: 



d 

+Uab {\b2j^i\^ + I62J+1P) a2j, (1) 

^ A t 2 

*^^^2j + l = —Jo^ * + {a2j + 0.2j+2) + £b,2j+lb2j + l + Ubb\b2j+l\ &2j + l + 



+ Uab (|a2jP + |a2j+2p) b2j+l, (2) 



where Aw = (wi— aj2), the parameters Uaa, Ubb define interaction between atoms 
induced by on-site atomic collisions and interaction between atoms of neighbor 
sites is defined by Uab- 

The first term in these equations takes into account the nearest- neighbor 
hopping induced by the stimulated Raman scattering. We assume that in- 
homogeneous broadening is absent, i.e., £a.2j = £a and £6,2^+1 = £&• If we 
introduce the control electromagnetic field amplitudes £1.2 then the nearest- 
neighbor hopping term read as Jq = ^12^1^2 j where fii2 is the matrix element 
of the Raman transition. If one introduce Jq — \Ja \ exp(iz?), and assume that 
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control electromagnetic fields have a constant phase, then i? can be included 
into complex value of the probability amplitudes 62^+1 • Thus we can substitute 
^2j+i ^ ^2^+1 exp(it9) — 62^+1, after that suppose the parameter Jq as real 
value. If the interaction between atoms of neighbor sites is neglected then the 
system of resulting equations takes the following form 

idh2j+i/dT = - {a2j +a2j+2) + 5b2]+i+ fib\b2j+i\^b2j+i, 

where Ae = {sb - Sa) ~ hAuj, (3^ = C/aa/| Jo|, = Ubb/\Jo\„ 5 = Ae/\Jo\. We 
use the normalized time variable t = t\Jo\/h. For the sake of simplicity we will 
assume that the exact resonance condition is hold, i.e., S — 0. 



3 Stationary solution 

It should remark that the atomic transport between neighbor sites in superlat- 
tice is absent if the phases of amplitudes a2j and 02^+2 as well as fo2j-i and 
^2j"+i will be opposite one. 

Fig. 1 represents schematically the probability amplitudes configurations 
in superlattice. Thin line arrows correspond to even sites the twin-line arrows 
correspond to odd sites. Sing of the probability amplitude is indicated by orien- 
tation of the arrow, i.e., plus (minus) corresponds to directed up (down) arrow. 
The configurations shown in Fig.l (a) and (b) are characterized by same energy, 
hence we can except the existence of the solution of the equations ([3]) which de- 
scribe the domain wall separated these two configurations. However, there we 
will not consider this case. 

The insertion of the ansatz a2j — (— l)-'a(r) and fo2j+i — i^^yb{T) into 
equations ([3]) results in following system of equations for the probability ampli- 
tudes 

^^=(3M'a, ^^^m'b (4) 
Equations ^ show that populations of the sites of each sublattice are indepen- 
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dent. It is convenient rewrite the equations Q in term of real variables 



a(T) = u(t) exp{v3„(T)}, 5(r) = w(r) exp{ip^{T)} 
that leads to the system of simple real equations 

Solutions of these equations read as 

u{t) = Uo, w{t) = Wo, faM = 'faO ~ Pa^O^ ^bi.T) ^ 'PbO ~ Pb^l (5) 

Choosing of the initial phases we can state the configuration of the initial prob- 
ability amplitude distribution as it shown in Fig.l (a)(<y9^Q = ipi^Q — 0), or in 
Fig- 1(b) {ip^o = 0,</Jbo = 7^)- 

4 Stability analysis for stationary distribution 

Stability of the solution found above will be analyzed in the framework of the 
linear stability theory. Let us consider the small perturbations of the stationary 
distribution 

= (-l)^a(T) + Sa2j, h^+i = (-1)^'6(t) + <562i+i, (6) 

with 

a(r) = Uo exp{-^^UQT}, 5(t) = wq exp{-/3fcWQr} 

The initial phases are chosen in the following form: Lp^^ = 0, <y9(,o = 0. 
The linear equations associated with ([3]) read as 
d 

d 

If one substitute (Sa2j(r) = p2j (''') exp{-/3^Mor}, (562j+i(t) = g2i+i (t) cxp{-/3^Wor}, 
than (l7|) and ([8]) can be rewritten as 
d 

^^^'zj = - [q2j-i + g2i+i) exp{i(^fc - <^„)} + /3„ug (p2i +P2i) > 
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d 

'^92i+i = - iP2j +P23+2) Gxp{i{ip^ ~ (^f,)} + Pb^l {q2j+l + q2j+l) ■ 
Assume that the constant probabiUty amphtudes (or the site population of 
superlattice) are related by the following expression 

/3a"o - Pb^l = Ai- (9) 

In this case the phase difference — ip^^ will be constant. We can put it to zero. 
Thus, the system of linear equations for small perturbations takes the form 

idp2j/dT = - {q2j-l + q2j + l) + Ai {p2j + P*2.j) , 

idq2j+i/dT = - {p2j + P2J+2) + Ai {q2j+l + q2j+l) ■ 
Substitution of the following expressions 

P2j = A exp{2ij kl) + B eyip{—2ij kl) 
q2j+i = Cexp{i(2j + l)fc/} + i:)exp{-i(2j + l)fcO 
into the differential-difference equations (flOl) leads to the system of linear dif- 
ferential equations 

idA/dr = -2 cos klC + Xi {A + B*) , 

idB/dr = -2coaklD + \i{A* +B), 

idC/dr = -2cosklA + Xi{D* +C), 

idD/dr = -2 cos klB + Xi{C* + D) , 
It is convenient introduce new variable ^ — Xit and constant parameter /i = 
2coskl/Xi = 2 cos fcZ//3^MQ. From the foregoing equations one can obtain the 
system of equations of second order 
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idA*/dT = 2 cos klC* - Xi {A* + B) , 

idB*/dT = 2cosfc;i:>* - Ai (A + B*) , 

idC*/dT = 2 cos kl A* -Xi{C* +D), 

idD*/dT = 2 cos klB* - Xi{D* +C) . 



-p^C + 2piA, 



(12) 



d^A/d^^ = 

a^B/af = -fi^B + 2fiD, d^D/d^^ = -fi^D + 2fiB 

Now, the characteristic equation for this system of equations (fT^ can be deter 
mined easily 

2/i ^ 

2/1 cr2 „ ^2 

_ ^2 2/1 

2/1 



T(ct) = Det 



v 



= 



(13) 
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Stability of the solutions of the equations ^ is determined by the roots of this 
equation, which can be written as 

4=/i2±2|/i| = (|^|±l)'-l. 

Instability of the configuration of site population under consideration means 
that imaginary part of the any root is not zero. But if (Tj_ is positive one, then 
ImCT — 0. One should note that o"^ > for any whereas cr^ > only at 
> 2. Hence, one can conclude that the configuration of site population is 
stable under following condition 

\ cos kl\>P^ul (14) 
In terms of physical meaning variable this inequality is read as 

Ubbwl = Uaaul < |Jo||cOSfcZ|. (15) 

Else, the modulation instability takes place if 

\ COS kl\ < (3 ^ul (16) 

The amplitude of small perturbations varies as exp(icr±^) ~ exp(i(T±AiT). As 
it was indicated above exponential growing of the amplitude is related with 
parameter ct_Ai. The imaginary part of a_Ai is the instability increment G{k), 
i.e., : 

G^{k) = 4|cosfc/| {P^ul -\coskl\) (17) 

If we consider the first Brillouin zone —tt/2 < kl < 7r/2, then the stability 
region lies into interval — arccos(/3^UQ) < kl < aiccos{f3 ^Uq) . The instabil- 
ity regions are determined by the inequalities — 7r/2 < kl < — arccos(/3^Uo), 
arccos(/3£jU§) < kl < tt/2. One can found that the instability increment is 
zero at boundary points of these regions. Maximum of the increment placed at 
points krnl = ± arccos(/3£(MQ/2) and maximum magnitude of increment is equal 
to Gm — G{km) = Pa'^o- ^^'^ P a'^Q > 2 the increment, ones taken at point 
kml = 0, has maximum value Gm — 2^ /3„Mq — 1. 
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Fig. 2 show the dependence of the instabihty increment on wave number of 
the weak perturbation and on nonhnearity parameter Ai = /3„Uq = Pi^Wq. This 
parameter is defined by the population of the superlattice sites. As one can 
see increasing of the Ai results in decreasing of the stability region. At Ai > 1 
the instability region occupies first Brillouin zone totally. All solutions (O are 
instable. 

5 Numerical analysis 

The aim of numerical simulation is study of the evolution of stationary con- 
figuration of the site population found above in response to a weak harmonic 
perturbation. Throughout this simulation relation I3^Uq = Pi^Wq is assumed. 
The equations ^ with 6 = were solved at the following initial (r = 0) condi- 
tions 



There we consider the periodic boundary condition: 6„+i = 6i, a„+2 — 0,2, 
where n is total number of sites in superlattice. The superlattice length was 
chosen to be multiple of the half-period of perturbation. 

As example we represent results of the numerical simulation of the modu- 
lational instability in superlattice containing 400 sites, where the initial values 
for probability amplitudes are a(0) = 2, 6(0) = 3, perturbation amplitude is 
6a = 0.01 and perturbation wave number is fc = 0.039/Z {I is distance between 
neighbor sites). If nonhnearity parameter Ai = /3(jWq is over one the modula- 
tional instability manifests itself causing an exponential growth of small per- 
turbations of the harmonic wave (Fig. 3) and (Fig. 4). By using the initial value 
probability amplitude of population for site of a -type, (i.e., a(0)) and the same 
value at time r, i.e. a(T), one can calculate the instability increment according 
to formula Gnum = lii[(a(T) — a(0))/(5a]/r (the same we can done for site of b 
-type). On the other hand, value of the increment G is determined by formula 
(dH). For (3^ul = 2.5 (that corresponds to = 0.625 and /3„ = 0.277) we 



a2j(0) 



(— l)-'uo exp{— iAir} + 5a cos{kl{2j)) 
{—lywo exp{— iAir} + 6acos{kl{2j + 1)) 
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obtained Gnum = 2 and G = 2.4. If the value /J^Mq is descried up to one under 
condition that all parameters of system are fixed, instability persists. However, 
the instability increment reduces progressively downstream. For P^Uq < 1 the 
regions of stability and instability appeared. Thus we can conclude that insta- 
bility state of the site population is typical for high populations of the sites in 
superlattice and for the case, where an on-site interaction between bosons dom- 
inates over photo-induced transport. Otherwise one can observe stable state 
picture: oscillation of an excess population near stable value, as shown in Fig. 5. 

Now we consider dynamics of the perturbations with different spatial fre- 
quencies. The nonlinearity parameter P^u^ , which put to be less than one, 
is fixed, but wave number of the harmonic perturbation k (see expression for 
initial conditions) will be varied. We can expect transition from stability to 
instability at kl = arccos(/3QMg) because of the stability region lies into in- 
terval - arccos(/3^Uo) < kl < arccos(/3aM§). Put /3^ul = 0.99 = 0.2475, 
/3f, = 0.11). Fig. 6 and Fig. 7 represent the time dependences of the site pop- 
ulations (i.e., square of modulus of probability amplitudes) respectively for 
fci/ = 0.10244 (cos(fci/) = 0.9948) and fciZ = 0.18124 (cos(/ciO = 0.9836). 

6 Conclusion 

In this work we have studied Bose-Einstein condensates in the one dimensional 
optical superlattice with two kinds of microscopic potentials (sites of the lat- 
tice). It was assumed that the deep of these potentials is enough the system to 
be described by the Bose-Habbard model. Steady state distribution of the prob- 
ability amplitudes of site population was found on base of the earlier derived 
system of equations ([3]) [22] , determining the dynamic of probability amplitude 
of sites filling (or the probability to find the boson in this site) in superlattice. 
Feature of this stationary state is the phase alternating, i.e., phases of prob- 
ability amplitude of one-type sites change on tt while going from one site to 
another. The stability of such state of BEG for small perturbations depends on 
the nonlinearity parameter /3qMq and also on the frequency of modulation by 
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itself. When the value of f3^UQ exceeds one, the found distribution is unstable 
for all wave numbers from Brillouin zone. However, when P^Uq is less than 
one the instability region in the Brillouin zone are determined by inequalities 
—tt/2 < kl < — aiccos{/3 ^Uq) and arccos(/3^UQ) < kl < 7r/2. So, the bosons dis- 
tribution is modulationaly instable for short wavelength perturbation, whereas 
the long wavelength perturbations dump out. If the nonlinearity parameter (or 
the average number of bosons in one site of superlattice) decreases, the instabil- 
ity region reduces. The values of the instability increment found by analytically 
were compared with the results of numerical simulation, and good agreement 
between them was found in the field, where the linear analysis of stability is 
valid. 

The nonlinear regime of the modulational instability was studied by using 
the numerical simulation. It was shown that the number of maxima in popula- 
tion distribution on sites appears. It should be remarked that the preliminary 
calculations demonstrate the development of chaotic behavior of the considered 
system. 

It is necessary to notice that in ordinary optical lattice the equations of mo- 
tion for the probability amplitudes of the sites population in continual limit may 
be transformed into the nonlinear Schrodinger equation having soliton solutions. 
It describes approximately the spatial solitons in ordinary optical lattices. In 
case under considering here the equations of motions ^ in continual limit result 
in more complex equations, which are not like to be completely integrable. So it 
is not necessary to expect formation of solitons chain as result of modulational 
instability. Inelastic interaction between solitary spatial waves appearing there, 
likely will leads to chaotic bosons distribution per sites. 

Finally, we assume it not required an especially effort for generalisation of the 
present results on the 2D cases, e.g., for the simple square or cubic superlattice. 
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FIGURE CAPTIONS 



Fig. 1. Two allowed configurations of the probability amplitudes distribution in 
superlatticc (a) Lp^^ = = 0, (b) Lp^^ = 0, (^,,0 = 
Fig. 2. Instability incrcnicnt G(k) versus k and Ai. 

Fig. 3. Strong rnodulational instability. Evolution of the site populations in 
superlattice for nonlinearity parameter I3^u^ = 2.5. . 

Fig. 4. Evolution of the sites populations for nonlinearity parameter Pau^ = 
2.5. Solid line corresponds to site with maximum of value of initial perturba- 
tion, dash line corresponds to site with maximum of negative value of initial 
perturbation, dot lines correspond to sites with intermediate value of the initial 
perturbations. 

Fig. 5. Same as in Fig. 4 but for nonlinearity parameter (3auQ = 0.8. and the 
dot lines correspond to sites with minimum value of initial perturbation. 

Fig. 6. Dynamic of weak harmonic perturbation with k\l = 0.10244, of the 
stability region.. 

Fig. 7. Same as in previous figure, but for harmonic perturbation with 
kil = 0.18124. One can see, that the distribution is instable now {Gnum = 0.08, 
G = 0.14). 
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